Wall-vortex composite solitons in two-component Bose-Einstein condensates 
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We study composite solitons, consisting of domain walls and vortex lines attaching to the walls in 
two-component Bose-Einstein condensates. When the total density of two components is homoge- 
neous, the system can be mapped to the 0(3) nonlinear sigma model for the pseudospin representing 
the two-component order parameter and the analytical solutions of the composite solitons can be 
obtained. Based on the analytical solutions, we discuss the detailed structure of the composite 
solitons in two-component condensates by employing the generalized nonlinear sigma model, where 
all degrees of freedom of the original Gross-Pitaevskii theory are active. The density inhomogeneity 
results in reduction of the domain wall tension from that in the sigma model limit. We find that 
the domain wall pulled by a vortex is logarithmically bent as a membrane pulled by a pin, and it 
bends more flexibly than not only the domain wall in the sigma model but also the expectation 
from the reduced tension. Finally, we study the composite soliton structure for actual experimental 
situations with trapped immiscible condensates under rotation through numerical simulations of the 
coupled Gross-Pitaevskii equations. 

PACS numbers: 03.75.Lm, 03.75.Mn, 05.30.Jp, 67.85.Fg 



I. INTRODUCTION 

Topological defects or topological solitons are solu- 
tions of systems obeying partial differential equations, 
representing localized structures with their stability be- 
ing due to non-trivial topology pj- Vortices in superfiu- 
ids/superconductors are an example of line topological 
defects 0, and it is believed that the analogous defects 
would exist in early universe as cosmic strings [3j. A do- 
main wall is a planer topological defect separating two 
different vacua or phases. When a symmetry group G of 
a system is spontaneously broken to a subgroup H, topo- 
logically allowed defect type is determined by the homo- 
topy properties of the order parameter space (vacuum 
manifold) G/H. In a (d + l)-dimensional spacetime, p- 
dimensional defects (p < d) exist if the homotopy group 
^d—p—i(G/H) is nontrivial. Thus, for d = 3 there will 
be planar defects (domain walls) if n (G/H) ^ 0, lin- 
ear defects (vortices or strings) if ■k 1 {G/H) ^ 0, and 
point defects (monopoles) if ^(G/H) ^ 0. These de- 
fects can be classified as "singular" or "continuous" in 
a sense whether (a part of) G is recovered at the core 
of defects or not. Order parameter is not defined at the 
core of a singular defect, while it is defined everywhere 
for continuous texture (defects). 

Bose-Einstein condensates (BECs) of ultra-cold atomic 
gases provide an ideal system for examining topological 
solitons in a quantum condensed system [J] . A major ad- 
vantage of this system is that the properties of BECs can 
be quantitatively described using the mean-field theory, 
namely, the Gross-Pitaevskii (GP) model. From experi- 



mental point of view, cold atom BECs are a versatile sys- 
tem to study topological defects, because most of the sys- 
tem parameters are tunable and optical techniques allow 
one to engineer the condensate wave function as well as 
to visualize the condensates directly. In the context of a 
single-component BEC characterized by scalar order pa- 
rameter with the broken U(l) symmetry, there are many 
papers discussing the properties of solitons and vortices; 
see Refs. [j| Q for reviews. In addition, realization of 
multicomponent (spinor) BECs with multiple order pa- 
rameters provides a ground to study more complex topo- 
logical solitons 0], as studied in superffuid 3 He Q. For 
example, the dark-bright solitons can be excited in two- 
component BECs [t| [l(| , where a dark soliton (density 
dip) of one component can trap a bright soliton (den- 
sity hump) of the other component . Exotic vortices 
composed of several order parameter components were 
observed experimentally |F2Hl6l ]. Because the order pa- 
rameter space of the multicomponent BECs can possesses 
higher symmetry than U(l) of the scalar BEC, their ho- 
motopy groups 7r„ with different n can become simulta- 
neously nontrivial and thus different kinds of topological 
solitons can coexist. There have been discussed the struc- 
ture, stability, and creation/detection schemes for various 
kinds topological solitons in multicomponent BECs, such 
as monopole [r7Fr2lj |. three-dimensional (3D) skyrmion 
[H-123, cosmic vortons (H 29], and knots 130, 113. 



In this paper, we discuss a 3D composite soliton con- 
sisting of domain walls and vortices in immiscible two- 
component BECs, as sketched in Fig. [TJ Two-component 
BECs have been realized by using the mixture of atoms 
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FIG. 1: (Color online) Schematic illustration of the wall- 
vortex soliton configuration in two-component BECs. The 
two-component BECs $1 at z > and 4>2 at z < are 
separated by the domain wall at the z = plane. Since 
penetration of the condensate densities takes place around 
the boundary, the domain wall is well-defined as a plane on 
which both components have the same amplitude, as shown 
in the middle panel. We assume that vortex lines are straight 
and perpendicular to the wall. Rotation is applied to the two 
components with the different rotation frequency fi, = fi^z 



with two hyperfine states of 87 Rb |32h35| or the mix- 
ture of two different species of atoms such as 87 Rb- 41 K 
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merits |35l . l38j demonstrated that miscibility and immis- 
cibility of two-component BECs can be controlled by tun- 
ing the atom-atom interaction via Feshbach resonances. 
The domain wall is referred to as a boundary of phase- 
separated two-component BECs and is well-defined as a 
plane on which both components have the same ampli- 
tudes (4p| - [43j . The vortices can be arranged by applying 
rotation of the confining potential around the z-axis to 
the phase-separated BECs [6fl . We assume that two com- 
ponents undergo phase separation in the z > and z < 
region, forming a domain wall lying at the z = plane 
and edges of the vortex lines along the z-axis attach to 
the domain wall. 

In our previous paper (44j, we pointed out that the 
wall-vortex composite soliton in two-component BECs 
can be identified as a non-relativistic analog of "Direchlet 
(DVbrane soliton" found in some field theoretical models 
[4a - l49| . This statement is based on the fact that the CP 
equations for two-component BECs can be mapped to 
the 0(3) nonlinear sigma model (NLcrM) by introdu cing 
a pseudospin representation of the order parameter [5CJ- 
[52j . The NLcrM admits the solitonic object that can have 
similar properties to the D-brane in the string theory pTjj . 
The purpose of this paper is to discuss in more detail 
the structure of this composite soliton in two-component 
BECs. The generalized NLcrM for two-component BECs 
includes additional degrees of freedom compared with the 
original 0(3) NLcrM, which modifies some properties of 
the composite soliton known in the previous literatures: 
(i) Vortices consisting of the composite soliton have sin- 
gular core, while they are nonsingular in the NLcrM. (ii) 



The density inhomogeneity of the BECs results in reduc- 
tion of domain wall tension from that in the NLcrM. (iii) 
The domain wall attached by a vortex is logarithmically 
bent, as a membrane pulled by a pin, and it bends more 
flexibly than not only the domain wall in the NLcrM but 
also the expectation from the reduced tension. We also 
study the composite soliton structure of rotating immis- 
cible BECs in a trapping potential through 3D numerical 
simulations of the coupled GP equations. To reduce the 
gradient energy of the density, the domain wall tends to 
be parallel to the rotation axis and forms a vortex sheet 
[EH . At high rotation frequency, a lattice of 2D skyrmion 
can form upon the domain wall, which undergoes trian- 
gular or rectangular ordering caused by the effective in- 
tercomponent repulsion realized in the restricted system 
on the domain wall. 

This paper is organized as follows. In Sec. [TH we formu- 
late the problem for two-component BECs and introduce 
the pseudospin representation to reduce the GP model 
into the NLcrM. In Sec. IIII1 we examine the structure of 
wall- vortex composite solitons based on the analysis of 
the NLcrM, where analytic solutions of these solitons can 
be obtained. In Sec. IIV1 we discuss how the composite 
solitons in two-component BECs are modified from the 
analytic solutions in the NLcrM and presents the results 
of 3D numerical simulations for the trapped immiscible 
two-component BECs under rotation. We conclude this 
paper in Sec. [V] 



II. THEORETICAL FORMULATION OF 
TWO-COMPONENT BECS 

We study the detailed properties of the composite 
solitons in two-component BECs, whose basic config- 
uration is illustrated schematically in Fig. [1] Two- 
component BECs are represented by the order parame- 
ters (*!, * 2 ) T = (y/Pie idl , y/p^e W2 ) T , which are the con- 
densate wave functions with the density Pj and the phase 
Oj (j = 1,2). They are confined in some trapping poten- 
tials and undergo phase separation, which results in the 
domain walls. The quantized vortices can exist in each 
component, being created by rotating the system or im- 
printing the circulating phase by the atom-laser coupling 
[fjj. We first show that the theoretical formulation of 
this system can be mapped to the NLcrM. This mapping 
was firstly discussed in the two-component Ginzburg- 
Landau theory for charged two-component Bose systems 
54( | . which was applied to the two-component BECs by 
some of the authors [EH [EH ■ 

The solitonic structure in two-component BECs is 
given by the analysis of the two-component GP model. 
The energy functional is given by 
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Here, rrij and pj are the mass and the chemical potential 
of the jth component, respectively. The trapping poten- 
tial is written by an axisymmetric harmonic oscillator as 



(2) 



with an aspect ratio A, where A < 1 (> 1) represents 
a cigar-shaped (pancake-shaped) potential. The coeffi- 
cients <7u, (722, and g±2 represent the atom-atom interac- 
tions. They are expressed in terms of the s-wave scat- 
tering lengths an and 022 between atoms in the same 
component and <zi2 between atoms in the different com- 
ponents as 



9jk 



2"kTi CLjk 
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where we have introduced the effective velocity field 
v c ff = v g + v s — 2A coming from the gradient of the 
total phase: 



h 



2m 



(7) 



and the flux flow of the spinor: 



v. = 7T- £(c;vo-Gvc;) 

Imi A — ' J 

h s z 
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Here, we have also used the relation 



with mj k = rrij + m k . The vector potential A is 

generated by (i) the rotation of the system A = (ft x r)/2 
Q or (ii) a synthesis of the artificial magnetic field by 
the laser-induced Raman cou plin g between the internal 
hyperfine states of the atoms [621 ] . 

The two-component GP model can be transformed to 
the similar form of the NLcrM by introducing the pseu- 
dospin representation of the order parameter. Here, we 
confine ourselves to the simple situation with the equal 
mass 771 1 = 777,2 = m and equal trapping frequency 
tjj\ = L02 = u)] its derivation in the case of the general 
parameters of the system, e.g., the mass imbalance and 
the difference of the trapping frequencies, was considered 
by Mason and Aftalion [52|. The condensate wave func- 
tions are denoted as 
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Here 
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C = [Ci,C2] T is the spin- 1/2 spinor with |(i| 2 + 
IC2I 2 = 1. The four degrees of freedom of the original 
wave functions = ^/pJe lBj (their amplitudes pj and 
phases 6j) are expressed in terms of the total density 
p = pi + P2, the total phase = 9\ + 6* 2 , and the polar 
angle 9 and azimuthal angle <f> of the local pseudospin 
s = (s x ,s y ,s z ) defined as 



, (5) 



where cr is the Pauli matrix, cos# = (pi — p2)/p, 4> = 
62 — Oi, and |s| 2 = 1. By using these variables, the to- 
tal energy Eq. ([T]) can be rewritten as the form of the 
generalized NLcrM [Sj: 



s = CV = 
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The second term in the right hand side of Eq. © corre- 
sponds to the classical NLcrM for Heisenberg ferromag- 
net. The generalized NLcrM has several unique features 
that are revealed as : (i) There is a gradient term of the 
total density, (ii) The spin stiffness, a prefactor of the 
(Vs a ) 2 term, is dependent on the total density p and 
is generally spatially inhomogeneous. (iii) There is an 
additional kinetic-energy term mpv^ s /2, associated with 
the presence of the superfiuid velocity v e g 7^ and the 
external vector potential A^O. 

The potential U is a function of the total density p 
and the ^-component s z of the pseudospin only, being 
explicitly written as 



U{p, s z 



c + CiS 2 



c 2 s 2 



(10) 



with 



Co 

Cl 
C2 



~[p(Sn + 522 + 2g 12 ) - 4(mi + Ma)], (H) 



P 



[p{m -922) - 2{pi - p 2 ) 



yCs 1 " +522 - 2512). 



(12) 
(13) 



If gn 7^ g22 7^ 912 or p\ 7^ P2, the anisotropic terms 
with the coefficients c\ and C2 break the global SU(2)- 
invariance of the system. The coefficient Ci can be in- 
terpreted as a longitudinal magnetic field that likes to 
align the pseudospin along the z-axis. The term with 
the coefficient C2 determines the spin-spin interaction as- 
sociated with s z ; it is antiferromagnetic for C2 > and 
ferromagnetic for C2 < [51(. The stationary point of 
this potential gives the equilibrium values 



(.922 - 9\2)P\ + (.9ii - .912)^2 



(.92 



911S22 
■512)^1 - 



Sl2 



(.9ii + 512V2 



(.922 - .912V1 + (.911 - .912)^2 ' 



(14) 
(15) 
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The determinant of the Hessian at that point is given by 



Im u 



d 2 Ud 2 U ( d 2 U 
dp 2 ds 2 \ds z dp 



[(.922 - 312W + (gn - 912)^2]' 

4(311322 - 3l2) 



The stationary point is a minimum or a maximum only 
when 311322 — 3i2 > 0. Otherwise, the minimum of 
the potential disappears within the range — 1 < s z < 1 
and the degenerate energy minima are given by s z — 1 
or s 2 = — 1. This situation corresponds to the ferro- 
magnetization, namely, the phase separation of the two- 
component BECs, which is discussed in the following. 



Reu 




III. TOPOLOGICAL SOLITONS IN 
NONLINEAR SIGMA MODEL 

In order to understand the properties of wall-vortex 
(D-brane) soliton in field theoretical model, we review 
the work by Gauntlett et al. [45] . The NLerM is a scalar 
field theory whose (multi-component) scalar field defines 
a map from a 'space-time' to a Riemann (target) mani- 
fold. The massive hyper-Kahler sigma model employed 
by Gauntlett et al. corresponds to the massive NLcrM 
for the effective description of the Heisenberg ferromag- 
net with spin-orbit coupling. The energy functional is 
given as 



E[s] = i 



dr 



£(V Sq ) 2 + [/(s) 



,0=1 



(17) 



also known as the Landau-Lifshitz model governing the 
high-spin and long wavelength limit of ferromagnetic ma- 
terials. Here, the amplitude of the vector is |s(r)| = 1 ev- 
erywhere. The ground state is two-fold degenerate such 
as s z = +1 and —1, where the potential is now described 
as 



U( S )=ml(l-sl 



(18) 



with a mass parameter m a . 

Under several conditions, our model Eq. © can be re- 
duced to the same form of Eq. (jTTJ) 44] . For the simple 
situation, we consider the homogeneous system without 
the trapping potential V = and put the parameters as 
311 = 322 = 3 and p\ = jj,2 = p. The anisotropy coef- 
ficient ci in Eq. (ITU1) then vanishes and the total energy 
can be written as 



E 



dr 



2m 



(VVP) 2 



mp 2 
-7T v ctt 



9 M 
-2V--9 



912 2 



p'(l-si) ,(19) 



where the constant term has been omitted. The coeffi- 
cient of the last term in Eq. (TT9"]) is positive because we 
consider the case of 312 > 3, giving the mass for the s z - 
field. For the limit 3 — > 00, which corresponds to the 



FIG. 2: (Color online) Stereographic projection from the 
sphere to the tangent plane at the north pole. 

Thomas- Fermi limit Q , we can approximate that the to- 
tal density is frozen to be p — fi/g = po and the (Vy/p) 2 - 
term vanishes. The kinetic energy associated with the 
superflow v c g is assumed to be negligible for simplicity 
|55j . By using the healing length £ = h/^/2mgp as the 
length scale, the total energy reduces to 



E = 



E 



gpie ' 

with the mass 



(Vs«) 2 +m 2 (l- S 2 ) 



1 



912 

9 



(20) 



(21) 



Therefore, the following discussion based on Eq. (fl7|) 
can be applied approximately to our system. Actually, 
as seen later, the additional degrees of freedom in two- 
component BEC system yield only quantitative modifi- 
cation of the soliton structure. 

To this end, we introduce the stereographic coordinate 



Sx T" 2 t S^ 

l + s z 



(22) 



where u = (00) corresponds to the north (south) pole 
of the spin sphere, as shown in Fig. [21 Then, each com- 
ponent of the pseudospin is written as 



(S X , Sy , S Z ) 



and Eq. ([TT]) becomes 



E= d 



(1 



\2\2 



(23) 



(24) 



The solutions of the topological solitons can be gained 
by taking the Bogomornyi-Prasad-Sommerfield (BPS) 
bound for the total energy (5(| Often, by insist- 

ing that the bound is satisfied (called "saturated"), one 
can obtain a simpler set of partial differential equations 
to solve, the Bogomol'nyi equations, from a square root 
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completion. Solutions saturating the bound are called 
BPS states and their energy is proportional to a topo- 
logical charge that characterizes the solitons. Here we 
summarize the properties of the BPS saturated solutions 
of the topological solitons. 



A. Vortex 

In the case of m a — in Eq. (IT7)) . the hamiltonian of 
the system has 0(3) symmetry. Since the symmetry of 
the ground state is broken to 0(2), the order parame- 
ter space is G/H = 0(3)/0(2) ~ S 2 . Then, the second 
homotopy group is nontrivial as 7r 2 (5 2 ) = Z. This sug- 
gests the presence of point-like defects such as monopoles 
and two-dimensional non-singular defects such as "2D 
skyrmions" (coreless vortex) [15|, [l6[ , because the former 
configuration can be mapped to the latter through the 
stereographic projection. 

First, we derive the analytic solutions of the coreless 
vortices by taking the BPS bound. We restrict ourselves 
to consider static solutions which are translationally in- 
variant along the z-axis. The total energy can be written 
as 



2 d x ud x u* 



dfx 



OyUdyU* 



= 2 



(i + M 2 ) 2 

\d w u\ 2 + \d w u\ 2 



(1 + M 2 ) 2 



(25) 



where we have introduced w = x + iy, d w = {d x — id y )/2 
and dw — (d x + id y )/2. The topological charge T v is 
given by the topological degree of the map u: R 2 — » S 2 . 
By considering the normalized area element of S 2 , the 
degree of u is given by 



2^ 



du A du* 



1 



52 (1 + M 2 ) 2 7T 



\d w u\ 2 - \d w u\ 
' {l + \u\ 2 ) 2 



•(26) 



The topological charge of the 2D skyrmion is given by 
T v 6 Z = iT2{S 2 ), giving the winding number of vortices 
passing through a certain z = const, plane. Inserting 
Eq. (j2"6")l to Eq. (|2"5) , we find that the total energy can 
be written as the sum of the topological charge and a 
positive correction: 



E = 2nT v 



d 2 x- 



\0a 



1 



|2\2 ' 



(27) 



Thus the energy is bounded by the topological charge 
E > 2irT v and the equality holds only if 



d,j,u = 0. 



(28) 



This equation is called a Bogomol'nyi equation. It is a 
first order equation whose solution gives a field configu- 
ration with a minimal energy within a fixed topological 
sector T v . Equation (|28l) also shows that u is a holomor- 
phic function of w only. Note that u is allowed to have a 



4- . s . 

2- 



till''' 

III/'- 
\ \\ 1 I / / ' > 

- VN ^ it"" 



— 1 v^.-»-» — 

\ \ \ \ v - 



■// / 




FIG. 3: (Color online) The spin profile of a coreless vor- 
tex. The spin component is given by (s x , s y , s z ) — (2x/(l + 
r 2 ),2y/(l + r 2 ),(l-r 2 )/(l + r 2 )). 



pole at any point w = Wi because its image on the tar- 
get S 2 is just the north or south pole. The requirement 
that the total energy is finite, together with the bound- 
ary condition that u has a definite limit as \w\ — > oo, 
forces u to be a rational map: 



u v (w) 



fH _ n 

g{w) n 



(29) 



where / and g are polynomials in w with no common 
factors. This solution gives the vortex configuration, in 
which f(w) and g(w) represent N n vortices (north poles) 
and N s antivortices (south poles), respectively. The po- 
sitions of the vortices are denoted by w™ and wf. Note 
that the total energy does not depend on the form of 
the solution, but only on the topological charges. In the 
NLerM, the energy is independent of the vortex positions 
w™' s ; in other words, there works no static interaction 
between vortices. 

Figure [3] represents the profile of the spin field for the 
simple vortex solution u v (w) — w. The spin orients up- 
wards at the center (s z — > 1 as — > 0) and it con- 
tinuously rotates from up to down as it moves outward 
radially (s z — > —1 as |iu| — > oo). The spin configura- 
tion of this continuous (coreless) vortex is known as a 
lump in field theory [58| , an Anderson- Toulouse vortex 
in superfluid 3 He [59(, or a 2D skyrmion in spinor BECs 
15, 161. 



B. Domain wall 

Next, let us consider the case for m a ^ 0. Assume that 
the field configuration s interpolates between two degen- 
erate ground state, so that s z = —1 for z — > — oo and 
s z = 1 for z — > oo. Then, we can have the domain wall 
between the two ground states. For such a configuration, 
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we obtain the BPS bound of the energy as 
d z ud z u* + ml\u\ 2 



(a) 



E 



dz 
dz 



(i + M 2 ) 2 

\d z u + m a u\ 2 m a (u*d z u + ud z u*) , 



(1 + M 2 ) 2 



(i + M 2 ) 



2V2 



;so) 



Here, the second term corresponds to the charge of the 
domain wall as 



T 



+oc (u*d z u + ud z u*) 

-oo Z (1 + |«| 2 ) 2 

— / dzd z s z = — [sj^^ 



(31) 



Under the above boundary condition, the charge becomes 
T w = +m CT . The energy is bounded by the domain wall 
charge E > T w and the saturated solution satisfies the 
Bogomol'nyi equation 

d z u + m a u = 0. (32) 

Then, we can obtain the BPS wall (kink) solutions 

u w = e -mA*-*o)-i<h j (33) 
or, in terms of s z , we have 

s z = tanh[rn CT (z — zq)]. (34) 

We can also consider the related solution with T w — 
— m ff , called anti-wall (antikink), and is obtained by mak- 
ing the replacement u — > vT x [s z — > —s z ). Here, zq repre- 
sents the position of the flat domain wall (s z — 0) whose 
transverse shift causes the Nambu-Goldstone mode due 
to breaking of the translational invariance. The phase 
—4>o corresponds to the azimuthal angle of the pseu- 
dospin s, causing the breaking of the global U(l) sym- 
metry locally along the wall. By promoting these two 
variables to dynamical fields as zq — > Zo(x,y,t) and 
4>o — > <Po(x,y,t), we can construct an effective theory 
of the domain wall. In the relativistic context, the low- 
energy dynamics of a single domain in the NLcrM wall 
can be described by the DBI action [lij |4(| , where the 
local U(l) gauge fields living on the wall are created by 
the dual transformation of the localized zero mode of the 
phase <pQ. This is the important ground why the domain 
wall in the NLcrM can be identified as an analog of a 
D-brane [11111 El. 



C. Wall-vortex complexes: D-brane solitons 

By combining the above two solutions of the topologi- 
cal solitons, we can construct the solutions in which vor- 
tices and domain walls are coexist. For a fixed topological 
sector, namely, for vortices (a domain wall) parallel (per- 
pendicular) to the 2-axis, the total energy Eq. ([24]) can 
be bounded with the topological charge as 
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FIG. 4: (Color online) The simplest wall- vortex soliton, where 
a single vortex along the 2-axis for z < is connected to the 
domain wall. The spin profile is explicitly given by Eq. (|38[) 
with m CT = 1. (a) The profile of s z . The spin texture of (b) 
z — plane and (c) y = plane. In (c), the wall position 
given by Eq. (f39)l is shown by solid curve. 



The BPS states can be represented with the form of the 
separated coordinate variables [13] 



it(r) = u v (w)u w (z) 



ntr(- 



w n ) 



,(z-zo)-i4 



(36) 



where u v and u w are satisfied with the Bogomol'nyi equa- 
tions (|28|) and (|32|) . and the forms of the solutions are 
given by Eq. (gSJ) and Eq. ([53]). 

To see the properties of the typical solutions, we depict 
the profile for the solutions of a single vortex and a single 
wall, written as 



u(w, z) 



(37) 



Here, we choose z = 0, 0o = 0, N n — 1, N s = 0, 
and to™ = 0. This is the simplest composite wall-vortex 
solution. The corresponding spin profile is written as 




\1 





1 






2ye~ m ° z 


1 


F \w\ 2 e- 2m " z 


1 


- \w \ 2 e ~ 2m ° z 


1 + \w\ 2 e- 2m " z 



\ 



J 



(38) 



E > IT, 



■ 27r|T v 



(35) 



which is shown in Fig. |4l For fixed |to| = r we have a 
domain wall solution along z-direction but for fixed z we 
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The wall position is given by 



-4 -2 2 4 



-4 -2 2 4 



FIG. 5: (Color online) The D-brane soliton to which two vor- 
tices attach, (a) The isosurface of s 2 = for the solution 
Eq. (|40p of the NLctM. The corresponding spin textures s in 
the z = plane and y = plane are shown in (b) and (c), re- 
spectively. In (c) , the wall position given by Eq. (|4ip is shown 
by solid curve. The wall becomes asymptotically flat due to 
the balance between the tensions of the attached vortices. 



have a vortex configuration in the x-y plain. For fixed 
z the coreless vortex has a scale size exp(m cr z); thus the 
size becomes infinity (zero) as z — > +oo (z — s- — oo). The 
wall position, i.e. the isosurface of s z = 0, is described 
by a logarithmic function 



1 



■ In liol 



(39) 



as seen in Fig. 0Jc). This situation is equivalent to the 
logarithmic bending when a membrane with a tension 
T is pulled by a pin at r — 0, where the profile of the 
membrane is given by z = (In r) /T as a problem of me- 
chanics. In the BPS solution, the tension of the domain 
wall is T w = m a as shown above. 

We can construct solutions in which an arbitrary num- 
ber of vortices are connected to a domain wall from 
Eq. (|36|) . because of the absence of the static interaction 
between vortices. Figure [5] shows a solution in which two 
straight vortices along the z-axis are connected to the 
wall on both sides. We assume that the vortices have 
ends at the positions w = xq and — xq upon the wall. 
Then, 



u{w, z) 



Xq 



Xq 



(40) 



z = In 



x 



xo 



(41) 



which becomes asymptotically flat (z —> 0) for |io| — > oo. 

It is instructive here to understand the wall configura- 
tion Eqs. (|39|l and (14"T1) in terms of tension of the vortices 
and walls schematically. The wall bending can be inter- 
preted to be caused by the tension of vortex attached to 
the wall. For the case of Eq. (|4Tj) , the tensional forces 
by the two vortices balance and the equilibrium position 
of the wall is well-defined as z = asymptotically. How- 
ever, for Eq. (j3"9")l the position is ill-defined since the force 
is unbalanced in the presence of a vortex only in the one 
side of the wall. Similarly, the wall has an equilibrium 
position z = zq when the number of vortices in one side 
equals to that in the other side, namely N n — N s in Eq. 
(f36f . The stiffness of the wall is represented by the coef- 



ficient Tn~ in Eqs. (|39|) and (|4T~|) with the wall tension 
T w = m a . Therefore, the wall is more flexible as the wall 
tension m a decreases. 

It should be mentioned here why this wall- vortex com- 
posite soliton has been referred to as "D-brane" soliton in 
the relativistic theory [45|, |46[ . On a single D-brane, the 
Abelian gauge theory is realized. The present domain 
wall has a localized U(l) Nambu-Goldstone mode and 
it can be rewritten as the U(l) gauge field on the wall, 
which is a necessary degree of freedom for the Dirac- 
Born-Infeld (DBI) action of a D-brane. Gauntlett et al. 
|45j have shown that Eq. ([55)1 reproduces the "Blon" so- 
lutions of the DBI action for D-branes in string theory 
by constructing an effective theory of the domain-wall 
world volume with collective coordinates Zo(x,y, t) and 
<fio(x,y,t) in u w (z), where </>o is periodically identified as 
4>a — > <Pq + 2tt. In the relativistic theory, the low energy 
effective action for these collective coordinates is given 
by 



I = -T K 



^yJ-detiGij+di^djM, (42) 



where (£°, £\ £ 2 ) = {t,x,y), and = %j + d t z^djZ 
with rjij — diag(— 1, 1,1) is the metric induced from the 
Minkowskii metric for a deformed membrane. Using the 
localized phase <f>o, we can introduce the U(l) gauge field 
Aj by taking a dual as 



di<j>o = eijkdjAk- 



(43) 



The effective action of zq and At corresponds to the so- 
called DBI action of the D2-brane: 



-T 

1 w 



cP£\/-det(G 4 , + F t 



(44) 



where Fij = diAj 



djAi is the electromagnetic field 
strength. The solution of this effective theory in the 
background of a point source with an electric charge and 
a scalar charge is known as Blon and its profile is pre- 
cisely coincident with that of the wall-vortex soliton in 
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the NLcrM (45|. Thus, the endpoints of the vortex lines 
in the NLcrM can be seen as electrically charged particles 
within this effective theory EH , and the domain wall 
can be seen as a D-brane on which fundamental strings 
terminate. However, the correspondence should be mod- 
ified in our non-relativistic theory which should be con- 
sidered in more detail but is beyond the scope of this 
paper. 

IV. WALL- VORTEX COMPOSITE SOLITON IN 
TWO-COMPONENT BECS 



the global shape of the condensate is elongated along the 
rotation (z) axis, because such a configuration minimizes 
the interface area between the two domains to decrease 
the energy cost due to the surface tension. We thus pre- 
pare a cigar-shaped trap with A = 1/4 to reduce the 
interface area and to keep the interface parallel in the x- 
y plane. We fix the intra-species s-wave scattering length 
as a = 5.61 nm and consider that the inter-species one 
a\2 is a free parameter in the following. The use of the in- 
terspecies Feshbach resonance will be crucial for realizing 
such a situation experimentally [Hj]. 



The mapping into the NLcrM can allow one to identify 
the domain wall of the two-component BECs as a non- 
relativistic counterpart of the D-brane soliton. Based 
on the analytic solutions of topological solitons in the 
simplified NLcrM, we next consider the structure of the 
wall- vortex composite soliton in trapped two-component 
BECs. The generalized NLcrM Eq. © has additional 
terms which are absent in the original NLcrM. Here, 
we discuss the modification of the soliton structure in 
the two-component BEC from the analytical solutions of 
Eq. ®. 

For simplicity, we assume the symmetric parameters 
mi = rri2 = m and g\ = g2 = g = 4iTh 2 a/m. We in- 
troduce the external trapping potential of Eq. ©. To 
nucleate and stabilize the vortices in trapped conden- 
sates, the system is supposed to be rotated at a rotation 
frequency Jl = f2z. In order to compare the numerical 
results with the previous analytic al results directly, we 
introduce the length scale £ = %/ ^/2mgp(0) and the en- 
ergy scale p = gp(0), where p(0) is the total density at 
the center of the trapping potential and can be estimated 
easily by applying the Thomas- Fermi approximation 
The coupled GP equation derived from the energy func- 
tional of Eq. ([T]) in the rotating frame of can be written 
as 



-V a +V , + |* i | 2 +7|**| 2 -afii, 



(45) 



where a — (Ait p(Q)aa\ Q ) 1 with cih = yh/muj and 
7 = c/ 12 /c/. The wave function has been scaled as 
^ — > J / o(0) , 5. In the following, we confine ourselves 
to the parameter range 7 > 1 for which the phase sepa- 
ration occurs. The trapping potential V can be written 
as 



2 V a ho, 



(46) 



and the rotation frequency is O = f2/wj_ The wave func- 
tions are normalized as J dr\^j\ 2 = Nj/[£ a p(0)]. The 
numerical solutions shown below are calculated by the 
imaginary time propagation of Eq. (jH} . 

Note that, in order to realize the configuration as 
shown in Fig. [I] the resulting domain wall should be per- 
pendicular to the rotation axis. Then, it is desirable that 



A. Domain wall 

We first discuss the structure of a domain wall in two- 
component BECs on the basis of the results of the NLcrM. 
The domain wall structures have been also studied by 
several authors [42|, EH . The position of the domain wall 
is defined as a plane on which two components have the 
same amplitude — \^f 2 \ ( s z — 0). Because the trap- 
ping potential does not play an essential role in the do- 
main wall structure, we consider the homogeneous sys- 
tem with V = as well as Ct — 0. Then the system is 
characterized by only one parameter 7 = g\ 2 /g- 

Let us assume that the wall lies in the z = plane and 
impose the following boundary conditions (^ r i,4 , 2) T — > 
(1,0) T at z -> +00, and (*i,* 2 ) T -> (0, 1) T at z -> 
—00. The typical profile of the domain wall solution is 
shown in Fig. [6] As z increases from negative to positive, 
the amplitude of $2 decreases as it approaches to the 
domain wall, while that of increases from zero with 
being apart from the domain wall. By representing the 
solution with the total density p = P1+P2 and the density 
difference, i.e. the ^-component of the pseudospin s z = 
(pi — p%)/ p, it can be clarified that the domain wall has a 
two-component structure, as shown in Fig.^b). The two 
length scales can be derived from the generalized NLcrM 
Eq. (|19j) . its dimensionless form being given by 



E = dr 



a 

1 



2 , ~~ 2 

pv e s 



+ -{p-l) 2 + -ml~p 2 (l-sl) 



(47) 



Here, we put p = p/p{0) and v e g = (mt;/h)v e ff; the tilde 
is omitted in the following discussion. By assuming that 
the system is uniform in x- and y-directions and v e g = 0, 
we consider the system spatially-dependent only on the 
^-direction. Using the identity 



-^-^[(Vs*) 2 + {s y Vs x - s x Vs y ) 2 }, (48) 

1 s „ 
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FIG. 6: (Color online) The numerical solutions of a domain 
wall in two-component BECs for 7 = 2 with the boundary 



conditions ^2 



(1, 0) at z +00, and (*i, * 2 



(0, 1) at z — > —00. (a) The density profile of pi (solid curve) 
and P2 (dashed curve), (b) The total density p — pi + P2 
(dashed curve) and the z-component of the pseudospin s z 
(solid curve). There are two characteristic length scales £ p 
and £ s in the solution. The results of the variational ansatz in 
Eqs. (|52[) and (|53[) for the profile of s z and p are also plotted by 
dotted curves, where A p in Eq. (|53[) is a variational parameter. 



we can write the energy of the problem as 

E z = dz 



2 2/1 2\ 



(49) 



The stationary solutions of the system satisfy the equa- 
tions 



ps z (d z s z ) 2 (d z p)(d z s z ) + p(d 2 s z ) 



(I -si) 2 



1 



sl) + P = 1,(50) 



m^ 2 s z = 0.(51) 



The asymptotic form of the profile can be obtained by lin- 
earizing with respect to p and s z around the ground state 
value p = 1 and s z = ±1 as p ~ 1 — e ±z ^ p = 1 — e ±v ^ z 
and s 2 ~ =Fl ± e ±z /£ s = =Fl ± e ±v/ ^ m " 2 for z — > =Foo. 
This gives the characteristic length scales £ p = 1/ v2 and 
^ s = l/(y2mo) in unit of £. Similar two-component 
structure can be seen also in the vortex solutions of two- 
component BEC [63| . In the strongly segregating limit 



FIG. 7: (Color online) (a) The relation between m'^ and 7. 
The circles are obtained by the fitting of the numerical solu- 
tion and Eq. (|52[) with the fitting parameter m^. We show 
the result for the cr-model limit m a — — 7I by the solid 
curve for comparison, (b) Tension of a domain wall as a func- 
tion of m a — ^/|1 — t|. The solid curve represents the results 
calculated from the ansatz of Eqs (|34[) and (|53p . where the 
energy is optimized with respect to the variational parameter 
A p . The dots show the results calculated with numerical so- 
lutions. The dashed line represents the tension T w = m CT in 
the cr-model limit 



g±2 — > 00, the domain wall is characterized by a single 
length scale £ p because £ s vanishes. 

Generally, the domain wall solution in the NLctM is 
written as Eq. (1551) . In terms of the condensate wave 
function, the domain wall solution can be written as ^1 = 
f ld (z - z )e-^°/ 2 and * 2 = f 2d (z - z )e^°/ 2 , where 
fjd is a real function with the wall center zq and the 
phase 0o, which can be identified as the relative phase 
between two components 4>q — <f> — 62 — 6\. The fixing 
of zq is due to the breaking of translational invariance by 
the given wall solution, while (f>o is due to the breaking 
of global U(l) around the domain wall, i.e., a narrow 
overlapping region of the two-component wave functions, 
and consequently there appears a U(l) Nambu-Goldstone 
mode localized around the wall. This feature satisfies a 
part of the requirement discussed in Sec lIII Bl and IIII CI 
that the domain wall in two-component BECs can be 
referred to as a D-brane soliton [44j |. 

It is instructive to consider the analytical form of the 
domain wall solution. The domain wall in the NLctM 
of Eq. pp|) under the g — > 00 (£ p —> 0) limit has a sin- 
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gle characteristic length £ s , where the profile is given by 
Eq. (1331) or When the spatial gradient of p is small 
enough for 7 > 1, the profile of s z in the generalized 
NLcrM must follow that in the NLcrM. We find that the 
domain wall solutions in Fig. [6] also follow correctly this 
profile function with slightly modified mass m' a : 

s z = tanh(m' z). (52) 

We make a fit of the numerical solution to Eq. (|52|) to 
extract the fitting value of m' a , which is plotted as a 
function of 7 in Fig. Efa). The mass parameter is al- 
most in agreement with the values given in Eq. (121[) in 
the sigma model limit, although it is slightly deviated as 
7 increases. Thus, the domain wall in two-component 
BECs can be regarded as the same solitonic object in the 
original NLcrM, including the quantitative details of the 
structure. 

The remaining total density can be described by the 
ansatz 



p = 1 — Apsech 




where A p is a variational parameter. According to Fig. 
Ef a )j h is reasonable to put m! a s» m a in Eq. ([52]) as 
s z = tanh(m CT z) = tanh(z/\/2^ s )j in other words, s z is 
assumed to be given by Eq. (|34l) with zq = 0. Inserting 
Eqs. ([M]) and (JS3]) into the energy of Eq. (gSJ) and mini- 
mizing the energy with respect to A p , we obtain the do- 
main wall profile semi-analytically as shown in Fig. [BJb) . 
The ansatz of Eqs. (f34|) and (|53l) agrees with the numer- 
ical result almost perfectly. The optimized energy corre- 
sponds to the tension of a domain wall E™ ln = T w , which 
is an extended version of Eq. (f3~Tj) for the two-component 
BECs. The tension is simply given by T w = m CT in the 
sigma model case, while it is significantly reduced in the 
BEC case because of the additional contribution of the 
total density, as shown in Fig.[7|b). 

B. Axisymmetric structure of wall- vortex complex 

We next consider the axisymmetric wall- vortex soliton 
in trapped two-component BECs. As shown in Fig. [TJ 
the (^2) domain is placed at z > (z < 0) and each 
component is assumed to have a straight vortex line at 
the center. The axisymmetric solution tyj — fj(r, z)e mj6 
with the real function fj, the polar angle 6 and the vor- 
tex winding number rij can be obtained by numerically 
solving the coupled CP equations: 

[ / d 2 Id d 2 nA ~ .a . 2 1 . 

= /i-(54) 

The parameters are the ratio of the coupling constants 
7 = 312/3 and the winding number rij. Here, we assume 



for simplicity that both components have the same par- 
ticle number N, setting m = ms7 Rb , ui — 20 x 27r Hz, and 
N = 10 5 . 

Before proceeding the discussion, we give some notes 
on the numerical solutions. In the energy-minimization 
process of the numerical simulations, the chemical po- 
tential is usually fixed in a homogeneous problem with- 
out a trapping potential, so that the particle number of 
each component is not conserved. Then, the pressure 
difference between two components, originated from the 
asymmetry of the solution, leads to the decrease in the 
population of the energetically unfavorable vortical com- 
ponent during the imaginary time evolution. Eventually, 
the vortex-free component fills all space as a final equi- 
librium solution. To obtain the desired solution, we have 
to adjust the chemical potential difference to balance the 
pressure between the two components, which is a trou- 
blesome task. Thus we make the numerical minimization 
by fixing the particle number in each component in the 
presence of the trapping potential, which is an experi- 
mentally relevant situation. 

Figure [5] (a) and (b) show the profile of the total den- 
sity p = p\ + P2 for 7 = 2 and the positions of domain 
wall, namely s z = (|Wi| = l^l) for several values of 
7, respectively. The vortex in the ^-component near 
the domain wall forms a coreless vortex, where its core 
is filled by the density of the ^1 -component and trans- 
forms into a singular vortex with distance from the do- 
main wall. Thus, the total density p is reduced at the 
position of the domain wall and vanishes at the singular 
vortex core around r = for z < —zq. The appearance of 
the singular core can be understood from the generalized 
NLcrM. When we give the phases 9\ = and 8-2 = for 
the case of Fig. OJa), for example, we obtain 

i 2 s = \{V9 2 ) 2 {1 - s z ) 2 = ^(1 - Sz ) 2 - (55) 

This kinetic energy density vanishes in the "fi {s z = 1) 
domain, while it contributes to the energy as pr~ 2 in the 
^2 (^ = ~ 1) domain. The latter divergent contribu- 
tion makes the singular vortex core in the ^-component 
around r = 0. 

This inhomogeneity of p implies that the assumption of 
the uniform total density to derive the NLcrM in Eq. (|20|) 
is not good and the solution is expected to be deviated 
from Eq. (|3"T]) . Nevertheless, the spin texture of this so- 
lution is almost identical to that in Fig. [4] The plot of 
the wall position in Fig. [5] (b) can be well fitted by the 
logarithmic function z = Ahir + B, as expected from 
Eq. (|3T)|) . Thus, the qualitative structural feature of the 
wall- vortex composite soliton in two-component BECs is 
not changed from the BPS solutions of the NLcrM. Ac- 
cording to Eq. (p9|) . we apply an analogy of a pulled 
membrane to this situation and extract the effective ten- 
sion T e ff of the domain wall from the numerical fitting 
with the relation 

z - B = ^-ln\w\, (56) 

Jeff 
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FIG. 8: (Color online) The numerical solutions of the ax- 
isymmetric wall-vortex soliton obtained by the GP equation, 
(a) The profile of the total density for (711,112) = (0, 1) and 
7 = 2 (b) Semi-log plot of the wall position for several values 
of 7 = giijg- Each plot can be fitted by the logarithmic func- 
tion as z = A\ar + B with A = 2.59, 1.91, and 1.62 for 7 = 2, 
4, and 8, respectively, where the fitting is made by using the 
data for r > 5 to avoid the contribution of the singular vortex 
core, (c) The effective tension determined from the relation 
T cB = A' 1 (see Eq. l[gB))) as a function of 7. The circles are 
obtained from the numerical fitting. The solid curve corre- 
sponds to the cr-model limit: T e ff = m CT = yfl — 7I, while 
the dashed curve is domain wall tension obtained by the vari- 
ational approach, shown in the solid curve in Fig. Ob), (d) 
The profile of the total density for (m, 712) = (1, 1). 



Here, the influence of the vortex associated with v e g- is 
neglected. As shown in Fig.[5Jc), the value of the effective 
tension T c g = A -1 is significantly reduced from not only 
T w = m a = yjl— t| of Eq. (|20|) in the cr-model limit 
but also T w of the BEC domain wall in Fig. [T^b). This 
means that the domain wall in this composite soliton can 
be more flexible than that in a single domain wall. This 
further reduction of the tension may be attributable to 
the following effects: (i) The rotational flow of a vor- 
tex causes the density inhomogeneity, where the density 




FIG. 9: (Color online) The numerical solutions of the axisym- 
metric wall-vortex soliton obtained by the GP equation for 
7 = 2. The profile of the total density for (a) (ni, n^) = (0, 2) 
and (b) (711,722) = (0,3) (c) The profile of p\ at r = as a 
function of z for 712 = 1 (solid curve), 2 (dashed curve), and 
3 (dotted curve). 

changes as P2 ~ r 2 / (2 + r 2 ) around the singular vortex 
core. The density difference between p\ and P2 (pi > pi 
near the wall) can enhance pressure from >Fi to ^ and 
lead to bend the wall more flexibly, (ii) The trapping po- 
tential also gives rise to the density inhomogeneity and 
the pressure balance can be modified radially. The latter 
is probably minor effect because Fig. Efb) shows that the 
wall is well-fitted logarithmically even for large r. Also, 
the effective tension is almost equal to that of the solu- 
tion calculated in the uniform system, although p\ 7^ pi 

Figure [5Jd) represents the profile of the total density 
for (711,712) = (1, 1). Because of the balance of the vortex 
tension, the domain wall becomes flat. This situation cor- 
responds to xq = in Eq. P0|) in the sigma model. There 
are two singular vortices which has infinitely thin distri- 
bution (5-functional form) of the vorticity; thus we have 
only the domain wall structure because the relative phase 
between the two components is uniform everywhere. 

Figure [9] represents the solution for (721,712) = (0,2) 
and (711,712) = (0,3). In these cases, the size of the 
vortex core extends radially, and the core is filled by 
the "Fi-component to be identified as the coreless vortex. 
According to the BPS solution of the NLctM u(w, z) — 
e -m a z w n 2 ^ t ne position of the domain wall is expected 
to become 

z -B~^-ln\w\. (57) 

Jeff 

The logarithmic fitting z = Alnr + B (for r > 5) of 
these solutions shows A =2.59, 4.44, and 6.20 for rc.2 =1, 
2, and 3, respectively. This is fairly agreement with the 
property of the NLcrM solution but its increase is lower 
than the expected linear dependance. This means that 
the tensile force for a domain wall pulled by two vortices 
is weaker than that of two BPS vortices. 
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Note that in the case with ri2 > 1 in Fig. [9l the 
total density does not vanish at the vortex core, as 
seen in Fig. [9] (c). Since the core size of the multiply 
quantized vortex becomes large with increasing n% like 
~ [r 2 /{2 + r 2 )]™ 2 , the density of the vortex- free compo- 
nent can enter the core easily. On the other hand, the 
vortex core for ri2 — 1 is apparently singular without den- 
sity. This indicates that there is a critical core size that 
allows the filling of the density inside the core. From dif- 
ferent points of view, there is a critical ratio of the chem- 
ical potential (/i 2 / ' P\) c that determines whether the vor- 
tex core can be filled with the other non-vortex compo- 
nent for a given 7 [64| • The 2D simulation shows that the 
optimized (711,712) = (0, 1) vortex state for fix = /i2 = 1 
and 7 = 2 is actually characterized by the empty vortex 
core. 



C. Non-axisymmetric structure: a wall with 
multiple vortices 

Next, we remove the axisymmetric condition and cal- 
culate the equilibrium state by the imaginary time prop- 
agation of Eq. ([4"5"|) in full 3D space from a suitably pre- 
pared initial configurations. To realize the final equilib- 
rium configuration as shown in Fig. [TJ we prepare the 
phase separated state in which ^1 ("I^-domains with 
some phase singularities (seeds of vortices) are located in 
the z > (< 0) region as the initial state of the calcula- 
tion. 

The panels of Fig. [TU] show the 3D distributions of the 
density difference \pi — cx \s z \ of the equilibrium state 
for several J7; this presentation is suitable to visualize the 
region of the vortex core and the domain wall (surface of 
Pi — Pi)- This configuration is energetically stable since 
it is obtained by imaginary time propagation. For Q = 
0.40 we obtain the (rii, U2) = (1, 1) state. Contrary to the 
axisymmetric structure of Fig. HJc), the end point of the 
vortices in each component is spontaneously displaced 
from the center, corresponding to xq 7^ in Eq. (|40"| . 
While the energy is independent of x$ in the BPS solution 
of the NLcrM, this displacement is due to the fact that the 
vorticity should be distributed broadly near the domain 
wall so as to reduce the associated kinetic energy as well 
as to reduce the gradient of p. Also, the vortex line is 
slightly bent due to the elongated trapping potential [65[ . 
Because our calculation uses the same rotation frequency 
il for both components, the number of the nucleating 
vortices should be the same for the both components [66[ . 

When the rotation is further increased, multiple vor- 
tices form a lattice in each component. Then, the domain 
wall begins to incline from the z = plane [Fig. flOT b)- 
(d)] and eventually becomes parallel to the rotation axis 
[Fig. [Pole)], even though the interface area (energy) in- 
creases. This is a vortex sheet structure [53[. The rea- 
son why this vortex sheet structure is preferred is due 
to the fact that the absorption of the vortices into the 



domain wall leads to the decrease in gradient energy of 
the singular vortex cores. This effect is also absent in 
the composite solitons of the NLcrM, which is free from 
the density gradient energy. Actually, when we consider 
the Thomas-Fermi limit, the gradient energy of the vor- 
tex core decreases, so that the structure such as Fig. Q] is 
expected to persist. The example is shown mFig.CEOlT), 
where the particle number is three times larger than that 
of Fig. HOT e). In this parameter setting, the domain wall 
is nearly parallel to the z = plane 

Another interesting property of this system at high ro- 
tation frequencies is that an ordering structure of many 
interface defects can emerge due to the complicated in- 
teraction effect. In each domain far from the domain 
wall, singular vortices form a Abrikosov triangular lat- 
tice. However, singular vortices become coreless vortices 
near the domain wall and a lattice of 2D skyrmions forms 
on the domain wall. Typical example is shown in Fig. 1111 
It is important to notice that p\ = P2 on the domain 
wall and a miscible state is effectively realized in this 
restricted 2D system of the immiscible condensates. In 
this effective 2D system, the intercomponent coupling t/12 
may be modified as gff , which determines the lattice 
structure of the 2D skyrmions. Note that a lattice of 2D 
skyrmions prefers a square lattice to a triangular lattice 
[67| . In the parameter setting in Fig. [TTJ the vortex end- 
points are shifted relative to each other on the domain 
wall to form a rectangular lattice of 2D skyrmions. With 
increasing c/12, we can see that the tendency to form a 
rectangular lattice from a triangular lattice becomes re- 
markable by comparing (a) and (b). These features are 
consistent with the phase diagram of a vortex lattice in 
2D miscible two-component BECs [50, |52j, |67( , and may 
originate from the static vortex- vortex interaction }68| . 
which is absent in the BPS solution. Our numerical solu- 
tions show that the domain wall in Fig. ITTT b') is inclined 
from the z = plane in order to elongate shorter sides 
of the rectangles. This suggests that the inclination is 
not accidental but caused by the energetic constraint to 
realize a square lattice. 



V. DISCUSSION AND CONCLUSION 

We have shown that a wall-vortex composite soliton, 
referred to as a D-brane soliton in field theoretical mod- 
els, can be realized as an energetically stable solitonic ob- 
ject in phase-separated rotating two-component BECs. 
Based on the NLcrM derived from the two-component 
GP model, we obtain the analytic solution of the topo- 
logical solitons, such as domain walls, vortices, and their 
complexes, by taking the BPS bound of the total en- 
ergy, which is a widely used technique in the field theory 
[1|. The topological solitons in trapped BECs are found 
to have the almost same character with the BPS satu- 
rated soliution in NLcrM. The inhomogeneity of the total 
density modifies the profile of the soliton quantitatively 
through the reduction of the domain wall tension. The 
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FIG. 10: The configuration of a wall connecting multiple vortices in trapped BECs. The panels show the profile of the density 
difference \p\ — pi\ oc \s z \ (isosurface of \pi — p2 1 = 0.0005) for f2 = (a) 0.4, (b) 0.5, (c) 0.6, (d) 0.7, and (e) 0.8. The parameter 
values used are A = 1/4, TV = 5.5 x 10 4 , and a 12 = 2a. The lower panels show the top view of the upper panels. In (f), we 
show the equilibrium structure for Q — 0.8, A = 1/4, and a\2 = 2a but the particle number N = 1.65 x 10 5 which is as three 
times as that in (e). 



domain wall pulled by a vortex is logarithmically bent as 
the BPS wall in the NLerM, but it bends more flexibly 
than expected by the tension of the BEC domain wall. 
The numerical analysis of full 3D simulations reveals that 
the complicated energetic constraint has an influence in 
determining the equilibrium configuration, such as the 
surface tension of the wall, the gradient energy of the 
density, and interactions between vortices and those be- 
tween interface defects. The last statement opens a prob- 
lem how to consider the properties of an effective 2D 
system realized in an interface of multicomponent con- 
densates, which can be affected by the extra dimensions 
(bulk regions). 

It should be noted, however, that there is one signif- 
icant difference of the wall- vortex composite soliton be- 
tween BECs and the NLerM. In the BECs described by 
the GP model, the total density vanishes in the singular 
vortex core for 712 = 1 because the density of nonrotating 
component does not enter into the vortex core, as seen 
in Fig. Htc). The coreless vortex near the domain wall 
shrinks to a singular vortex for a finite distance and thus 
we can identify a point connecting a singular vortex and 
a coreless vortex. This is in contrast to the case in the 
NLerM, where a coreless vortex extends to infinity along 
the thin vortex core, avoiding the singularity since s is 
well-defined everywhere. Hence, a connecting point is ab- 
sent, or more precisely, it should be positioned at infinity 



in this model. In a field theoretical model, the connect- 
ing point forms defect called "boojum" , which serves the 
negative binding energy of vortices and a wall and a half 
of the negative charge of a single monopole [13, |48[ . Boo- 
jums are known as point defects existing upon the surface 
of the ordered phase; the name was first introduced to 
physics by Mermin in the context of superfluid 3 He (69j . 
Boo j urns can exist in different physical systems, such as 
the interface separatin g A and B phases of superfluid 3 He 
H [ZOl , liquid crystals 7l[ , the Langmuir monolayers at 
air- water interfaces [72 1. multi-component BECs with a 
spatially tuned interspecies interaction [73|, [74| , and high 
density quark matter [75j . In the present model, boo- 
jums can be found at the end points of vortices on the 
domain wall, at which the vortices change their charac- 
ter from singular to coreless type. A suitable topological 
charge for boojums in two-component BECs can be de- 
rived by noting the analogy of the Abelian gauge theory 
[76| . A detailed study of the distribution of the boojum 
charge and the interactions between boojums remains as 
a future study. 

As pointed out in Ref. [HJ], the domain wall in two- 
component BECs are useful to simulate some analogue 
phenomena of the D-brane physics in a laboratory. One 
famous example is a nonequilibrium dynamics such as 
brane-antibrane annihilation, which was proposed for a 
possible explanation of inflationary universe in string the- 
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ory. In braneworld scenarios of cosmic inflation the anni- 
hilation may lead to defect production that could be di- 
rectly observed in atomic BECs; the experiment has been 
performed with superfluid 3 He A-B interfaces, but the de- 
tection of defects is difficult |77j. Recently, we proposed 
that domain-wall annihilation in two-component BECs 
actually demonstrates a brane-antibrane collision and a 
subsequent creation of cosmic strings, causing tachyon 
condensation accompanied by spontaneous Zo symme- 
try breaking in a two-dimensional subspace |78j |. Also 
we propose that, when strings are stretched between the 
brane and the antibrane, namely when the filling com- 
ponent has vortices perpendicular to the wall, "cosmic 
vortons" can emerge via the similar instability [29|. All 
of these phenomena can be monitored directly in experi- 
ments. We hope that our works open a new trend of the 
cold atom physics as "simulator of everything" . 



FIG. 11: (Color Online) The equilibrium structure and the 
spin texture along the (s z ) = plane for fl — 0.8, A = 1/4, 
N = 1.65 x 10 5 , and (a) 012 = 2a and (b) a\% = 4a; the 
parameters in (a) are the same with those of Fig. 11014 ) . Left 
panels show isosurface of the density difference \p\ — p?\ as in 
Fig. 1101 Right panels show the distribution of the pseudospin 
field s(r) in the plane (a) along the z = plane and (b) 
slightly declined from 2 = plane, indicated by the bold lines 
in the left panels. The color scale shows the magnitude of 
s z . The circles and squares mark the position of 9 2- and V&i- 
vortices, respectively, where we only mark in the y < region 
for clarity. 
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